Parity quantum numbers in the Density Matrix Renormalization Group 
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In strongly correlated systems, numerical algorithms taking parity quantum numbers into account 
are used not only for accelerating computation by reducing the Hilbert space but also for partic- 
ular manipulations such as the Level Spectroscopy (LS) method. By comparing energy difference 
between different parity quantum numbers, the LS method is a crucial technique used in identifying 
quantum critical points of Gaussian and Berezinsky-Kosterlitz-Thouless (BKT) type quantum phase 
transitions. These transitions that occur in many one-dimensional systems are usually difficult to 
study numerically. Although the LS method is an effective strategy to locate critical points, it has 
been lacked an algorithm that can manage large systems with parity quantum numbers. Here a new 
parity Density Matrix Renormalization Group (DMRG) algorithm is discussed. The LS method is 
the first time performed by DMRG in the S=2 XXZ spin chain with uniaxial anisotropy. Quantum 
critical points of BKT and Gaussian transitions can be located well. Thus, the LS method becomes 
a very powerful tool for BKT and Gaussian transitions. In addition, Oshikawa's conjecture in 1992 
on the presence of an intermediate phase in the present model is the first time supported by DMRG. 

PACS numbers: 02.70.-c, 75.10.Pq, 05.10.Cc 



I. INTRODUCTION 

In onc-dimcnsional (ID) quantum many-body sys- 
tems, Bcrczinsky-Kosterlitz-Thoulcss (BKT) 1 and Gaus- 
sian type quantum phase transitions are usually diffi- 
cult to be studied numerically. New methods for ac- 
curate determination on the BKT and Gaussian criti- 
cal points have being proposed. 2-7 However, for exam- 
ples, precise detections from entanglement entropy need 
to compute on very large size N > 10000 systems 6 , and 
detections from bipartite fluctuations need the prelim- 
inary knowledge of Luttinger parameters 7 . The Level 
Spectroscopy (LS) method 8,9 based on the sine-Gordon 
theory is an old method for the BKT and Gaussian 
type quantum phase transitions. By comparing differ- 
ent excitation energies, the LS method is able to de- 
tect the Gaussian and BKT critical points accurately. 
Ground state phase diagrams have been studied by the 
LS method in many models. 10-14 However, implementa- 
tion of the LS method usually needs precise energy eigen- 
values with parity quantum numbers p = ±1. Although 
parity or space inversion is usually a good quantum num- 
bers in condensed matter physics, the Exact Diagonalza- 
tion (ED) 15 was the only accurate numerical algorithm 
with parity quantum numbers. Therefore finite size effect 
usually interferes the LS method because of the limited 
sizes in ED. In order to avoid the finite size effect in some 
difficult cases 14,16 , a precise algorithm with large sizes 
and parity quantum numbers is now desired urgently. 17 

Density Matrix Renormalization Group (DMRG) 18-21 
is a powerful numerical method of strongly correlated 
systems in ID and two-dimensional (2D) lattices. Lat- 
tice size in DMRG usually achieves hundreds of sites 
with high accuracy. Recent study on spin- 1/2 system 
in Kagome lattice revealed that DMRG may be the most 
accurate numerical method in frustrated 2D systems. 22 



Symmetry also plays an important role in DMRG. U(l) 
symmetry, such as conservation of number of particles 
and z-component of total spin, is the most frequently 
employed. Non-Abelian symmetry 23 and other descrete 
symmetries 24,25 are often more complicated and rather 
difficult to implement. For the parity symmetry, DMRG 
faces intrinsic difficulties in the algorithm. The core 
problem of parity in DMRG is that the algorithm divides 
supcrblock into left and right parts. This division de- 
stroys the spacial inversion if the length of left and right 
blocks are not equal. When they are equal, for exam- 
ple the infinite- system DMRG, parity can be utilized. 24 
However, finite-system DMRG or so-called sweeping pro- 
cedure is usually needed for pushing ground state into a 
reliable precision. In this case parity is difficult to be 
utilized. 
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FIG. 1. (a) Parity DMRG scheme proposed by S0rensen 24 . 
The left and right blocks must be the same size. Parity quan- 
tum numbers are only available in the infinite-system DMRG 
within this scheme, m is the dimension of a block and d is the 
dimension of a site, (b) Ladder scheme for parity quantum 
numbers in the sweeping procedure. The left dash rectangle 
represents the system-block and the solid rectangle represents 
the enlarged-system-block. The dimension of one site (dash 
ellipse) becomes d 2 . The bold bond controls the boundary 
conditions. Thus accurate results within open (OBC), peri- 
odic (PBC) and twisted boundary conditions (TBC) can be 
easily obtained. 
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In this paper, parity quantum numbers p = ±1 in 
the sweeping procedure are discussed in an alternative 
ladder scheme. The trick is very simple and no prior 
knowledge of group theory required. Therefore the LS 
method overcomes the limitation of small sizes in ED 
and becomes a very powerful tool in detecting BKT and 
Gaussian transitions. In the state-of-the-art of DMRG, 
matrix-product state (MPS) 26 is the main actor as vari- 
ational algorithms. 27 However, the language of MPS is 
not used to describe the trick, but traditional DMRG in- 
stead, since it continuously shows the significance in the 
studies. 22,28 In the following section, the trick is discussed 
in S=l/2 XXZ model. In Sec. Ill A, the proposed parity 
DMRG is applied to perform the LS method in S=2 XXZ 
model with uniaxial anisotropy. BKT and Gaussian tran- 
sitions arc both located precisely. Thus in Sec. IIIB, the 
boundaries of Intermediate- D phase are determined. In 
other words, the presence of Intermediate- D phase which 
Oshikawa predicted 20 years ago 29 is now supported by 
DMRG. The conclusions are summarized in Sec. IV. 
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of rung i. Thus V sys is the parity operator which 
only exchanges legs of the enlarged- system-block and 
[Psys,V sys \ = [S z sys ,V sys \ = [p S y S ,S z sys ] = 0. The eigen- 
states of reduced density matrix have both quantum 
numbers of S* and V sys which can be used for the next 
iteration in sweeping procedure. The bases of a rung, 
dashed ellipse in Fig. 1(b), | r< ) € {|ft>, \U), |4t>, \U)}- 
Each |tj) is an eigen-state of S z , however, these bases 
are not eigen-states of the parity Vi- Since [Vi, S z ] = 0, 
one can easily find the basis set which are mutual 

eigen-states of V% and Sf . 
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II. PARITY DMRG 



Consider the S=l/2 XXZ model, 



N 1 

H = 2^ 9 + ^i+i) + ^SjSj +1 , 



(1) 



where A is the anisotropic parameter and N is total num- 
ber of spins. In developing numerical methods of strongly 
correlated systems, this model can be regarded as a ba- 
sic testing model 15,18 because the dimension of one site 
d = 2 and it can be solved by using Bethe ansatz. 30,31 
The model undergoes a BKT quantum phase transition 
at A c = 1. Ground state energy per site at this critical 
point in the thermodynamic limit is eo = — In 2 + j . Par- 
ity quantum numbers of finite size ground state and first 
excited state depend on N. When L = N/2 is even (odd), 
the parity of ground state is also even (odd) and the first 
excited state is odd (even) . These two states become the 
two-fold degenerate ground state in the thermodynamic 
limit. 

The reason that quantum numbers of Sf ot can be 
utilized is the fact that reduced density matrix of 
enlarged- system-block p sys and the z-component of total 
spin of enlarged-system-block S* s are compatible, i.e., 
[psys, Sg y8 ] = 0. Thus each eigen-state of reduced density 
matrix has a quantum number of S z ys , and these quan- 
tum numbers can be used for the next iteration. How- 
ever the parity operator V is a global operator which only 
acting on the whole chain, the supcrblock. It is impossi- 
ble that in the original scheme to have an operator V sys 
which only acting on the enlarged- system-block. 

In order to make it possible, the ladder scheme should 
be employed. As shown in Fig. 1(b), the chain is folded as 
a two-leg ladder. Now the global parity operator V can be 



These states are so-called singlet and triplet states 
with corresponding quantum numbers of parity pi = 
—1,1,1, and 1, respectively, and quantum numbers of 
z-component of spin sf = 0,0,1, and — 1, respectively. 
Eq. (2) implies the quantum number of enlarged-system- 
block pk = piPn , where pi is quantum number of system- 
block and p n is quantum number of a single rung. These 
indexes satisfy the relation k = (I — l)d 2 + n. 
The {|cri)} forms a transformation matrix Q. 
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The spin operators of the first and second site of rung i 
are = Q T (S Z ® l dxd )Q and Sf l = Q T {l dxd g> S Z )Q, 



and so do S^, S^, and S^i- ldxd is a d x d iden- 
tity matrix. The two spins of the first rung interact 
with each other in the ladder scheme. The Hamiltonian 
of the first rung Hi is a d 2 x d 2 matrix and it is for- 
tunately already diagonalized with eigenvalue ~ 2 ~ A for 
singlet state and eigenvalues j and j for triplet 

states, respectively. The Hamiltonian of the latest rung 



TABLE I. Energy per site of ID spin-1/2 XXZ model. p and 
pi are parity quantum numbers of the ground state and first 
excited state, respectively. The number of state kept in the 
largest size is up to m = 1500 and the truncation error is of 
the order 10~ 9 . Three sweeps are performed. PBC is used 
here, i.e., Hi — Hl = diag(^A, 2=* A A) in E q. (5). 



N 


A 


po 


E /N 


Pi 


Ei/N 


54 


1 




-0.44343001 


+ 


-0.44189435 


100 


1 


+ 


-0.44322957 




-0.44277665 


162 


1 




-0.44317856 


+ 


-0.44300474 



3 



E/N 



-0.439 
-0.44 
-0.441 
-0.442- 
-0.443, 
-0.4441 



X E/N 
+ E/N 
. -ln(2)+l/4 



A' 







5x10 



1x10" 



FIG. 2. The ground state and first excited state energy per 



site of spin-i Heisenberg model, 



A = 1 in Eq. (5). PBC 



is used and three sweeps are performed in parity DMRG. The 
dot indicates the exact ground state energy in the thermody- 
namic limit from Bethe ansatz. 



Hl, the bold bond shown in Fig. 1(b), can be flexible 
chosen as different boundary conditions. For example, 
H L =H\= diag(^^-, is the periodic bound- 

is the open boundary con- 
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ary condition (PBC), Hl - 



dition (OBC), and H L = diag(^, |, |) is the 

twisted boundary condition (TBC). Thus Eq. (1) can be 
rewritten as 
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where L = N/2. Once a DMRG programmer has above 
information, the DMRG with parity quantum numbers 
can be performed. The energies of ground state and first 
excited state with corresponding parity quantum num- 
bers are listed in Tab. I and plotted in Fig. 2 for looking 
at convergence of energy with respect to system size. 

Although parity reduces the dimension of superblock 
by a factor 2, the disadvantage of the ladder scheme is 
the dimension of one site becomes squared. In other 
words, 4 in S=l/2 chain, 9 in S=l chain or t-J model, 
16 in Hubbard model and 25 in S=2 chain. The overall 
computational difficulty increases, especially when d is 
very large. One can further use single center site 20 to 
reduce the dimension of superblock. As mentioned by 
White, the ladder scheme is a better configuration for 
PBC but it only improves convergence with the number 
of sweeps. 20 I emphasize the feasibility of parity DMRG 
in the ladder scheme. Thus in the case d 2 = 25, the 
S=2 XXZ model with uniaxial anisotropy is examined in 
Sec. Ill as practical examples. DMRG with a single cen- 
ter site 20 does not be employed here. Besides quantum 
numbers (Sf 0t and V), the only optimization in this work 
is the wave-function transformation 19 which proposed by 
White in 1996 becomes a standard optimization in per- 
forming DMRG. 



III. S=2 XXZ MODEL WITH UNIAXIAL 
ANISOTROPY 



S=2 spin chains are current interesting research top- 
ics, including topological phases and quantum phase 
transitions 14 ' 32-37 , magnetization process 38 , and cold 
atoms loaded into a ID optical lattice 39 . The S=2 XXZ 
model with uniaxial anisotropy is defined by 
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where A and D are the XXZ anisotropy parameter and 
uniaxial anisotropy parameter, respectively. This model 
exhibits a Haldane gap in the Haldane phase with a non- 
local string order. 40,41 Although the ground state phase 
diagram for the entire parameter space is unclear, it is ex- 
pected to have Haldane, Neel, Large- .D, Intermediate-/?, 
XY1, XY4, and Ferromagnetic phases. Quantum phase 
transitions with various universality classes are in this 
model. It is still a controversial issue 14 ' 29,32,42-44 that 
whether the Intermediate- D (ID) phase which was con- 
jectured by Oshikawa 29 twenty years ago exists in ground 
state phase diagram. Although it was concluded absence 
of ID phase by DMRG 42-44 , surprisingly by using ED to 
perform LS on very small sizes N < 12, Tonegawa et al. 
obtained the ID phase in a very narrow region. However, 
the ED estimations may not reveal results of thermody- 
namic limit because of finite size effect, especially when 
D > 2. Energy level crossing does not appear in small 
sizes at all. The region of the phase diagram that can not 
be determined by ED was thus determined by extrapo- 
lation of phase boundaries. 14 



A. Level Spectroscopy method 

Therefore the first application of parity DMRG is to 
perform the LS method by focusing on the BKT tran- 
sition from XY1 to Large- D phase as well as the Gaus- 
sian transitions from Large-D to ID phase and ID to 
Haldane phase. Following the LS method, these three 
excitation energies Eq (N, , + ; tbc) , Eq (N, , — ; tbc) , 
and E {N, 2, +; pbc) should be compared. Where 
E (N,M,p; tbc) and E (N, M, p; pbc) denote the low- 
est energy eigenvalues of TV spins in the subspace of z- 
component of total spin M and parity quantum number 
p within TBC and PBC, respectively. 

The BKT pseudo-critical points of XYl-Large-D sat- 
isfy the condition 14 



£0(^,0,+; tbc) -E (N,2,+; pbc) = 0, 



(7) 



and the Gaussian pseudo-critical points of Large-73-ID 
and ID-Haldane satisfy the condition 14 



E o (N,0,+;tbc) -E o (N,0, -; tbc) = 0. 



(8) 



These pseudo-critical points arc extrapolated to N — > 00 
by assuming the quadratic scaling function of A^~ 2 after 
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FIG. 3. The LS method for BKT transition in the model 
Eq. (6). A = is fixed. Energy differences in Eq. (7) with 
different sizes are shown in (a) and the extrapolation of critical 
point is shown in (b). N=8, 10, and 12 are computed by ED, 
and larger sizes are computed by DMRG. Only N > 16 are 
used for the extrapolation. The XYl-Large-D critical point 
from the extrapolation D c = 2.796917 ± 5 x 1CT 6 . 



collecting the finite size results. Fig. 3 demonstrates the 
LS method for the BKT transitions. A = is fixed in 
Eq. (6). The number of states kept is up to m = 800, 
and truncation error is about 10 -8 . BKT quantum phase 
transitions usually difficult to precisely locate numeri- 
cally. However, this critical point is well located by the 
LS method after using the parity DMRG. As shown in 
Fig. 3(b), the critical point D c = 2.796917 ± 5 x 10" 6 
is determined precisely. By comparing the accuracy of 
recent development of detecting BKT quantum phase 
transitions, 7 the LS method becomes a very powerful tool 
for the BKT quantum phase transitions. 

For the Gaussian transitions, Fig. 4 shows the LS re- 
sults with fixed D = 2.1. Due to the finite size ef- 
fect, the energy level crossing does not take place until 
TV > 10. Consequently, it is not sufficient for quadratic 
fitting with only the ED data of N = 10 and 12. 14 In 
fact, the ID-Haldane pseudo-critical point A 2 behaves 
non-monotonically with size. As shown in Fig. 4(b) and 
(c), only 16 < N < 28 are used for the extrapola- 
tion. The critical points in the thermodynamic limit are 
A c i = 2.51610 ± 3 x 10~ 5 and A c2 = 2.53303 ± 4 x 10~ 5 
for Large- D-1D and ID-Haldane, respectively. The num- 
ber of states kept is up to m = 1300 in the largest size, 
and the dimension is about 5 x 10 7 in the Lanczos diag- 
onalization of finding the lowest energy eigenstate. The 
order of truncation error is 10~ 8 in the worst case. In 
contrast, recent accurate determination of the Gaussian 
transition achieved the similar accuracy by considering 
the entanglement entropy up to twenty thousand sites. 6 
Therefore the LS method also becomes a very powerful 
tool for the Gaussian transitions. In addition, following 
the LS method, it implies the presence of ID phase in 
model Eq. (6). Thus Oshikawa's conjecture in 1992 29 is 
now supported by DMRG. However, the region of the 
ID phase in D = 2.1 is merely about 0.017. It may ex- 
plain why the ID phase hasn't been found from previous 
DMRG studies consequently. 42-44 
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FIG. 4. The LS method for Gaussian transition in model 
Eq. (6). D = 2.1 is fixed. Energy differences in Eq. (8) 
with different sizes are shown in (a) and the extrapolations of 
critical points are shown in (b) and (c). N < 12 are computed 
by ED, and larger sizes are computed by DMRG. Since energy 
level crossing takes place when N > 10, it is not sufficient for 
quadratic fitting with only the ED data. The dash lines in 
(b) and (c) are linear guild from ED. The solid lines are the 
extrapolations with only 16 < N < 28. The critical points 
for Large- D-1D and ID-Haldane are A c i = 2.51610 ± 3 x 10~ 5 
and A c2 = 2.53303 ± 4 x 10~ 5 , respectively. 



B. Intermediate- D phase boundary 

The phase diagram is focused on the region D > 2.1 
which contains some ambiguity in previous studies. A 
possible scenario is that the ID phase boundaries of ID- 
Haldane and Large-!?- ID merge at a point. This point 
(A c , D c ) ~ (2.64, 2.19) was estimated by extrapolation of 
phase boundaries from D < 2. 14 In this scenario, the Hal- 
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FIG. 5. D = 2.14. Although the pseudo-critical points satisfy 
the condition Eq. (8), there is no energy level crossing in the 
thermodynamic limit. In other words, there is no Gaussian 
transition. This strong finite size effect is not observed 14 until 
using the new parity DMRG method. 



dane and Large- D phases are suggested to be the same 
phase. 14,37 However, from previous analysis by ED, 14 it 
is known that larger value of D leads the level crossing 
to begin at larger size. Therefore, based on the presence 
of the ID phase, there are two other possible scenarios 
against the first one. One is that the merge point (A c , D c ) 
exactly locates at the boundary of Neel phase and be- 
comes a multi-critical point. Another possible scenario 
is that the ID phase extends and reaches the Neel phase. 
In the latter two scenarios, Haldane and Large- D phases 
are completely separated by the ID phase and definitely 
not the same phase, though recent studies suggested they 
could be the same phase. 14,37 

Departing from empiric of ED, as shown in Fig. 5 for 
D = 2.14, despite the pseudo-critical points satisfy the 
condition Eq. (8) in small sizes, there is no energy level 
crossing in the thermodynamic limit. In other words, 
according to the LS method, there is no Gaussian transi- 
tion, and the scenarios that the ID phase completely sep- 
arates Haldane and Large- D phases are excluded. The 
more precise phase diagram by means of the new parity 
DMRG is shown in Fig. 6. 

Finally the accuracy for the current research is care- 
fully examined. It is known that the often-used measure 
of the error in DMRG calculations, the truncation error, 
may not reveal the true relative error of energy in the 
infinite- system DMRG. 2i In the finite-system DMRG, 



TABLE II. The energies for N = 28 in Eq. (8) as m increases. 
A = 2.51 and D = 2.1. Three sweeps are performed. 
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FIG. 6. The phase boundary of Intermediate- D phase in the 
S=2 XXZ chain with uniaxial anisotropy. the phase boundary 
lines between ID and Haldane phases and between Large-D 
and ID phases merge at a point between 2.13 < D c < 2.14. 



the truncation error is a more reliable indication to the 
true error. The energies of A = 2.51 and D = 2.1 for 
N = 28 in Eq. (8) as the number of states kept m in- 
creases are listed in Tab. II. From these data, the abso- 
lute errors of energy in Fig. 4(a) are argued to be of order 
10 -6 , and the relative errors to be the same order with 
the truncation error. Therefore, the energies are suffi- 
ciently accurate for distinguishing the energy differences 
in the current studies. The phase diagram in Fig. 6 thus 
provides a solid numerical evidence for the presence of the 
ID phase in the S=2 XXZ chain with uniaxial anisotropy 
Eq. (6). 



IV. CONCLUSION 

The parity (space inversion) quantum numbers are of 
significance not only in computational advances but also 
in the LS method. By employing the ladder scheme, the 
global parity operator can be decomposed into the prod- 
uct form in Eq. (2). It makes DMRG be able to uti- 
lize the parity quantum numbers. The LS method is the 
first time performed by DMRG in the S=2 XXZ model 
with uniaxial anisotropy. The BKT critical point of XY1- 
Large-Z) as well as Gaussain critical points of Large-D-ID 
and ID-Haldane are determined very precisely. Thus the 
LS method is suggested to be the most powerful tool for 
detecting BKT and Gaussian transitions. 

The phase diagram of the S=2 XXZ model with uniax- 
ial anisotropy, Fig. 6, is investigated by focusing on the 
ID phase boundary where the ED results are strongly af- 
fected by the finite size effect. This work consistent with 
previous finding 14,37 thus provides a first DMRG support 
to Oshikawa's conjecture 29 in 1992. 

Since the proposed method in the case d 2 = 25 works 
well, it is also feasible for a large number of mod- 
els, including S=l/2 spin chains, S=l spin chains, t-J 
model, Hubbard model, topological interacting fermion 
systems 16 , and ladder, etc. 
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FIG. 7. The old parity schemes in the infinite-system 
DMRG 24 for (a) PBC and (b) TBC. The dash line is the 
boundary coupling. The left and right blocks must be the 
same length. The truncation error may not reveal the true 
relative error of energy. 24 It is difficult to evaluate how large 
number of states kept m is enough. 
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FIG. 8. Relative error in the ground state energy for N = 100 
sites spin-1/2 Heisenberg chain as a function of the number of 
states kept m. £ D mrg are computed by finite-system DMRG, 
and the exact energy i? B cthe are obtained from algebraic Bethe 
ansatz. 31 
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Appendix: Boundary conditions and parity DMRG 

In the Level Spectroscopy (LS) method 8,9 , Twisted 
Boundary Conditions (TBC) and parity quantum num- 
bers are usually required in the numerical algorithm. Al- 
though the old parity infinite-system DMRG may ob- 
tain accurate energy within Periodic Boundary Condi- 
tions (PBC) as long as enough large number of states 
are kept, the PBC scheme in Fig. 7(a) can not be used 
within TBC because the parity described in Rcf. 24 is 
not conserved anymore. Instead of Fig. 7(a), the TBC 
scheme in Fig. 7(b) should be used. However it may not 
be a practical implementation because the direct con- 
nection between two blocks slows down the algorithm 
dramatically. 18 Moreover, the truncation error in the 
infinite- system DMRG may not reveal the true relative 
error of energy. 24 It is difficult to evaluate how large num- 
ber of states kept m is enough. These considerations may 
be the reason why previous LS studies were always per- 
formed by Exact Diagonalization. 10-14 

The scheme in Fig. 7(a) was first proposed by White in 
his initial DMRG papers for PBC. 18 It is known that the 
traditional DMRG performs worse for PBC. In order to 
improve the PBC case, the ladder scheme was first pro- 
posed by Qin et al. in 1995. 45 Unfortunately, the com- 
parison did not publish because the ladder scheme does 
not have a distinct improvement. Fig. 8 shows the com- 
parison of convergence with m for finite-system DMRG. 
While the new parity DMRG within PBC has merely a 
little superiority to the traditional DMRG, it should be 
emphasized that LS method which based on sine-Gordon 
theory has a strong numerical tool now. 46 For the PBC 
case, recent development by Pippan, White and Evertz 
have a remarkable improvement with the Matrix Prod- 
uct State (MPS) algorithm. Readers with interests are 
referred to Ref. 47 and 48 for more detail information. 
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